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Abstract. We propose a scheme involving cold atoms trapped in optical lattices to 
observe different phenomena traditionally linked to quantum-optical systems. The 
basic idea consists of connecting the trapped atomic state to a non-trapped state 
through a Raman scheme. The coupling between these two types of atoms (trapped 
and free) turns out to be similar to that describing light-matter interaction within 
the rotating-wave approximation, the role of matter and photons being played by 
the trapped and free atoms, respectively. We explain in particular how to observe 
phenomena arising from the collective spontaneous emission of atomic and harmonic 
oscillator samples such as superradiance and directional emission. We also show how 
the same setup can simulate Bose-Hubbard Hamiltonians with extended hopping as 
well as Ising models with long-range interactions. We believe that this system can be 
realized with state of the art technology. 
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1. Introduction 

Cold atoms trapped in optical lattices have been proved to be very versatile quantum 
systems in which a large class of many-body condensed-matter Hamiltonians can be 
simulated (see [TJ [3] for extensive reviews on this subject). One of the first proposals 
in this direction was the article by Jacksch et al. |4j, where it was proved that the 
dynamics of cold atoms trapped in optical lattices is described by the Bose-Hubbard 
Hamiltonian provided that certain conditions are satisfied; shortly after, the superfluid- 
to-Mott insulator phase transition characteristic of this Hamiltonian was observed in 
the laboratory |5J. Since then, the broad tunability of the lattice parameters, and the 
increasing ability to trap different kind of particles (like bosonic and fermionic atoms 
with arbitrary spin or polar molecules), has allowed theoreticians to propose optical 
lattices as promising simulators for different types of generalized Bose-Hubbard and 
spin models which are in close relation to important condensed-matter phenomena 
PQEJE]. Recent experiments have shown that optical lattices can be used to address open 
problems in physics like, e.g., high-T c superconductivity [6], to study phenomena in low 
dimensions such as the Berezinskii-Kosterlitz-Thouless transition [7], or to implement 
quantum computation schemes [E]. 

In this article we keep digging into the capabilities of optical lattices as simulators, 
showing how they can also be a powerful tool for simulating quantum-optical 
phenomena. In particular, following our previous proposal |9J, we introduce various 
schemes in which superradiance-like phenomena can be observed. In addition, we will 
show that these same schemes realize Bose-Hubbard and spin models with extended 
interactions (in particular, models where hopping between sites and ferromagnetic 
interactions between spins are long range instead of nearest-neighbors). 

The concept of superradiance was introduced by Dicke in 1954 when studying the 
spontaneous emission of a collection of two-level atoms [10] (see [H] for a review). He 
showed that certain collective states where the excitations are distributed symmetrically 
over the whole sample have enhanced emission rates. Probably the most stunning 
example is the single-excitation symmetric state (now known as the symmetric Dicke 
state), which instead of decaying with the single-atom decay rate T , was shown to 
decay with NTq, N being the number of atoms. He also suggested that the emission 
rate of the state having all the atoms excited should be enhanced at the initial steps 
of the decay process, which was a most interesting prediction from the experimental 
point of view, as this state is in general easier to prepare. However, Dicke used a very 
simplified model in which all the atoms interact with a common radiation field within the 
dipolar approximation; almost 15 years later, and motivated by the new atom-inversion 
techniques, several authors showed that dipolar interactions impose a threshold value for 
the atom density, that is, for the number of interacting atoms, in order for superradiance 
to appear [121 1131 E] • 

It is worth noting that together with atomic ensembles, spontaneous emission of 
collections of harmonic oscillators were also studied at that time [13] . Two interesting 
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features of this system were reported: (i) an initial state with all the harmonic oscillators 
excited does not show rate enhancement, on the contrary, most excitations remain within 
the sample in the steady state, while (ii) the state having all the oscillators in the same 
coherent state has a superradiant rate. 

In the last years, there has been a renewed interest in this phenomenology because 
of its potential for quantum technologies [151 EH El HI HI 120] • In these new approaches, 
superradiant effects differ from that found by Dicke in that they are mediated by the 
initial entanglement of the emitters, not by a build up of their coherence. Of particular 
interest to our current work is the analysis performed in [18J (see also [1 TJ ) of the 
spontaneous emission by regular arrays of atoms separated a distance do comparable to 
the wavelength A of the emitted radiation. It was shown that directional emission can 
be obtained for particular single-excitation entangled states when A > 2do, the emission 
rate being enhanced by a factor x ^ (A/<io) 2 iV 1//3 . 

Also of interest for our purposes is the work carried by Sajeev John and collaborators 
on the collective emission of atoms embedded in photonic band-gap materials [2T1 122]. 
In this kind of materials the density of states of the electromagnetic field is zero for 
frequencies laying within the gap, what gives rise to the phenomenon of light localization 
[23J. As for collective emission, it was shown by using a simplified Dicke-like model 
[211 122] that when the atomic transition lies close to the gap, radiation is emitted 
in the form of evanescent waves and the single-excitation symmetric state has a rate 
proportional to iV 2 r , which is even larger than the one predicted by Dicke in free space. 

We will show that all this phenomenology can be observed in optical lattices with 
the setup already presented in [9], which is depicted in figured! Consider a collection 
of bosonic atoms with two relevant internal states labeled by a and b (which may 
correspond to hyperfine ground-state levels, see for example [1]). Atoms in state a 
are trapped by a deep optical lattice in which the localized wavefunction of traps at 
different lattice sites do not overlap (preventing hopping of atoms between sites), while 
atoms in state b are not affected by the lattice, and hence behave as free particles. A pair 
of lasers forming a Raman scheme drive the atoms from the trapped state to the free one 
[1], providing an effective interaction between the two types of particles. We consider 
the situation of having non- interacting bosons in the lattice [24], as well as hard-core 
bosons in the collisional blockade regime, where only one or zero atoms can be in a 
given lattice site [25]. In the first regime, the lattice consists of a collection of harmonic 
oscillators placed at the nodes of the lattice; in the second regime, two-level systems 
replace the harmonic oscillators, the two levels corresponding to the absence or presence 
of an atom in the lattice site. Therefore, it is apparent that this system is equivalent 
to a collection of independent emitters (harmonic oscillators or atoms) connected only 
through a common radiation field, the role of this radiation field being played by the free 
atoms. This system is therefore the cold-atom analog of the quantum-optical systems 
considered above, with the difference that the radiated particles are massive, and hence 
have a different dispersion relation than that of photons in vacuum. Moreover, we show 
that the field of free atoms can be characterized by a dispersion relation which is similar 




Figure 1. Scheme of our proposed setup. Atoms in state a are trapped in an optical 
lattice, while atoms in state b are free, and can thus have any momentum. An external 
pair of Raman lasers connect the two levels with some detuning A. We will show 
that above some critical value A = A c the atoms in state b are able to leave the trap 
(left), while below this, they are trapped forming a bound state with atoms in state 
a (right). This behavior is typical of band-gap systems where the density of states is 
zero between the connected states. 

to the one obtained for photons within a photonic band gap material [HI [2TJ [22] . It is 
then to be expected, and so we will prove, that this system will show the same kind of 
phenomenology as its quantum-optical counterparts. 

To conclude this introduction let us explain how the article is organized. We first 
introduce the model and find its associated Hamiltonian (Section [2]). Then we study 
the radiative properties of the system when the lattice sites emit independently (Section 
|3]); we will show in this section that our system shows localization of the free atoms, in 
analogy to the localization of light occurring in a photonic band gap material, and we 
identify as well the Markovian and non-Markovian regimes of the emission. In Section 
H]we study collective effects. We first deduce a reduced master equation for the lattice 
atoms, explaining under which conditions it is valid. We then show how by changing 
the system parameters extended Bose-Hubbard and spin models (Section I4.ip . Dicke 
superradiance of atomic (Section 14.2. ip and harmonic oscillator (Section 14 . 2 . 2 [) samples, 
and directional superradiance (Section I4.3P can be observed. We finally talk about the 
effect of restricting the motion of the free atoms to 2D or ID traps in Section [5j and 
give some conclusions in Section [6j 

2. The model 

Our starting point is the Hamiltonian of the system in second quantization [26] . We will 
denote by \a) and \b) the trapped and free atomic states, respectively (having internal 
energies froj® and fajjl). Two-body interactions for the trapped atoms are included with 
the usual contact-like pseudopotential [3j [24] , but we neglect the collisions for the free 
atoms. The Hamiltonian is then written as H = Hq + H a _b, with 
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H a - b = hnj dW^™^ (r) (r) + H.c; (lb) 

H contains the individual dynamics of the atoms, Hj being the first-quantized motion 
Hamiltonian of the atom in the corresponding state, and g = 47rH 2 a s /m, where a s is the 
s-wave scattering length of the trapped atoms (which have mass m). H a - b contains the 
Raman coupling between the atomic states, k L = k x — k 2 and ojl = oj\ — u 2 (laser wave 
vector and frequency in the following) being the relative wave vector and frequency of 
the two lasers involved in the Raman scheme (see figured]), with Q the corresponding 
two-photon Rabi frequency. 

For atoms in state \a), H a = — (h 2 /2m) V 2 + V opt (r), where V opt (r) corresponds 
to a 3-dimensional optical lattice with cubic geometry and lattice period do- We work 
with ultracold atoms under conditions such that their wavefunctions can be described 
by the set of first-band Wannier functions localized around the nodes of the lattice [I]. 
The traps of the optical lattice are approximated by isotropic harmonic potentials PQ, 
what allows us to write the Wannier functions as 

(r - rj) = exp [- (r - rj ) 2 /2X 2 ] , (2) 

where ij = d j is the position of the j 6 Z 3 lattice site (we consider M sites in each 
orthogonal direction defining the cubic lattice), and Xq = h/muo, being w the frequency 
of the harmonic trap. The energy associated to these wave functions is Eo = 3huo/2. 
We will assume that Wannier functions localized at different lattice sites do not overlap, 
hence preventing tunneling between sites. 

On the other hand, atoms in state \b) can move freely in every direction of space 
according to = — [h 2 /2m) V 2 , and hence the plane waves ijj^ (r) = e lk ' r /VV, with 
energy £ k = H 2 k 2 /2m, are their motion eigenfunctions (V is the total available volume 
for the free atoms, which we might take as infinite for calculations). 

We consider two opposite regimes for the interaction between trapped atoms 
(remember that in (11 qp any two-body interaction involving the free atoms has been 
neglected). The first limit consists in neglecting the interactions, which might be 
accomplished by, e.g., tuning the scattering length with an additional magnetic 
field through a Feshbach resonance [24J. Expanding the quantum fields onto their 
corresponding motion eigenstates, we get within the interaction picture: 

H a - b = Y,9k^ kt -^ L)r >^l + H.c, (3) 



with 



g k = Ml\l — exp 



ix 2 (k-k L ) 2 



^-A, (4) 
2m K ' 



A = ul — (u)b — w«) being the detuning of the laser frequency respect to the \a) ^ \b) 
transition (u a = w° + 3u /2 and uj b = The operators {aj,aj} and {b^, bl} satisfy 
canonical bosonic commutation relations, and create or annihilate an atom at lattice 
site j and a free atom with momentum k, respectively. 
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As for the second limit, we assume that the on-site repulsive atom-atom interaction 
is the dominant energy scale, and hence the trapped atoms behave as hard-core bosons 
in the collisional blockade regime, what prevents the presence of two atoms in the same 
lattice site [25] : this means that the spectrum of ajaj can be restricted to the first two 

states ||0)j , j, having or 1 atoms at site j, and then the boson operators jaj, aj j 

can be changed by spin-like ladder operators |<7j , Oj j = {|l)j(0|, |0)j(l|}. In this second 
limit the Hamiltonian reads 

H a - b = ^^e^'-^-^^Oj&t + H.c (5) 

Hamiltonians fl3]) and (p]) show explicitly how this system mimics the dynamics 
of collections of harmonic oscillators or atoms, respectively, interacting with a common 
radiation field. Note that, particularly, the dispersion relation appearing in (jl]) is similar 
to that of the radiation field within an anisotropic 3D photonic band-gap material, where 
photons acquire an effective mass close to the gap j2TJ |22]. We will show how by varying 
the parameters of these Hamiltonians, the system has access to regimes showing many 
different phenomena, as explained in the Introduction. 

Note finally that in order to satisfy that the trapped atoms are within the first 
Bloch band, it is required that u 3> A, Q. 



3. Emission of an atom from a single site 

Before studying the collective behavior of the atoms in the lattice, it is convenient to 
understand the different regimes of emission when the sites emit independently. To this 
aim we first study the emission properties of one atom in a single site following the 
analysis performed in [211 E2] for an atom embedded in a photonic band-gap material. 
We are going to show that there exists a critical value of the detuning above which the 
atom is emitted, while below which the radiated atom gets bound to the trapped atom, 
and there is a nonzero probability for the atom to remain in the trap (see figure [1]). 
This second regime is the analog of the photon-atom bound state predicted more than 
25 years ago [23]. We will also identify the Markovian and non-Markovian regimes of 
the emission. 

Following the Weisskopf-Wigner procedure [271 EH] , we write the state of an atom 
in a single site (which we take as j = 0) as 

|V (t)> = A (t) |1, {0}) + B * (*) 1°' M > (6) 

k 

where |1, {0}) refers to the state with an atom in the trap and no free atoms, and |0, Ik) 
to the state with no trapped atoms and one free atom with momentum k. We can use 
the Schrodinger equation to write the evolution equations of the coefficients A and B^; 
then, the equation of B^ can be formally integrated arriving to a single equation for A 




Figure 2. (a-c) Time evolution of the population when the lattice sites emit 
independently, (a) corresponds to the trapping region, (b) to the pure non-Markovian 
emission of the atoms, and (c) to the radiative region having both Markovian and non- 
Markovian contributions. The logarithmic scale reveals that the Markovian regime 
requires |A| ^> na 2 . (d) Steady state population for finite ojq as a function of the 
detuning A for different values of the Rabi frequency fi (see Appendix I). We can 
appreciate the phase transition from a trapped to a radiative state for A = 4Q 2 /u>o. 



given by 



A (t) = - / dt' G(t- t') A (t') 



(7) 



where we have assumed that the atom is initially in the trap, that is, A (0) = 1 and 
Bk (0) = 0. This is a non-Markovian equation where the free atoms enter the dynamics 
of the trapped atoms as a reservoir with correlation function 

^ £ \9k\ 2 e" iAfcT = T ^ ^ exp (.At) , (8) 



G\ 



3/2 



where we have taken the continuous limit for the momenta in the last equality and 
kk = for simplicity. Note that this correlation function coincides exactly with that of 
atoms in anisotropic 3D photonic band-gap materials [211 



In Appendix A we show how to handle equation (J7|) by using Laplace transform 



techniques. In particular, in Appendix A. 1| we show that in the strong confinement 



regime u — > oo, this equation can be solved as 

2a 



A{t) 



cc 



dx 



^— x + iA j + i47ra 2 x 



(9) 
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with a 2 = 8ft 4 /c^, A = A - 4fi 2 /w and 

c = 26+/ (6+ - &_) and 6 = 6+ for A < 0, 

c = for < A < vra 2 , (10) 

c = 26_/ (6_ - 6+) and 6 = 6_ for A > 7ra 2 , 



where 6+ = ^—1 ± y 1 — A/ira 2 

Solution ([9]) has a very suggestive form. The first term has exponential evolution 
as expected for a Markovian radiation process, while the integral term has a nontrivial 
time dependence which always goes to zero for sufficiently large times. The time while 
the integral term is comparable to the exponential one defines the non-Markovian part 
of the evolution. This result also shows that the behavior of the emission is highly 
dependent on the parameters of the system: 

• In the A < region Im {6 + } = 0, and hence the steady state population is nonzero 



and given by \A(t — > oo)| 2 = |c| 2 = f 1 — 1/yl + \A\/ira 2 j ; a fraction 1 — |c| 2 

is radiated during the non-Markovian period which lasts longer as A/a 2 goes to 
zero (see figure [2^,). In this region the radiated atoms are emitted in the form of 
evanescent waves, and the equivalent of a photon-atom bound state [23] is formed. 



In the < A < ira 2 region (figure [2b) the steady state population is zero, i.e., the 
atom leaves the trap eventually, and the evolution is dictated solely by the integral 
part of the solution (hence it is pure non-Markovian radiation). 

Finally, in the A > ira 2 region the atom is radiated in a Markovian fashion only for 



A ^> a 2 , with a decay rate given by T = Im {6 2 } ~ 4f2 2 y 2ttA/uq] if this is not 
the case (A ~ a 2 ), the integral part is comparable to the exponential part during 
most of the evolution time, and hence the radiation process is non-Markovian (figure 

no. 

Therefore, there exists a phase transition at A = 4Q 2 /uq in close analogy to that 
of spontaneous emission in photonic band-gap materials [21, 22J: Above this value the 
atom is forced to leave the trap, while below it there is a nonzero trapped population 
left in the steady state. Using the techniques we explain in Appendix A. 2 it is possible 
to show that this phase transition is also present for finite u (see figure [2H). 

Note that the Markovianity condition |A| S> a 2 can be recasted as |A| 3> T , which 
coincides with the usual Markov limit rg ^> tr [29J once the characteristic evolution 
times of the lattice atoms and the reservoir's correlation function are identified with 
rg = Tq 1 and tr = |A| _1 , respectively. 
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4. Collective dynamics 

Let us study now the collective emission of atoms from many sites. We analyze this 
through the master equation for the reduced density operator of the lattice within the 
Born-Markov approximation [29] . 

We introduce the Born approximation by fixing the state of the reservoir to a non- 
evolving vacuum. This leads to an effective interaction between different lattice sites 
mediated by the two-point correlation function 

G >-> W = h E^- 1 ^^^'^ 1 (ML) = -P ( ~f^ ~ ikx ■ G (r) , (11) 

being rj_i = rj — ri = do (j — 1). The Born approximation yields the dynamics of the 
system to the leading order in the coupling parameters [29] . On the other hand, as 
commented at the end of the previous section, the Markov approximation is valid only 
if |A| ^> r, where T is a typical evolution rate for trapped atoms, which can differ from 
To owed to a renormalization due to collective effects as we show below. We also include 
the single-site energy shift that we found in the previous section (the analogous of the 
Lamb shift in the optical case), and substitute A by A in the two-point correlation 
function. 

The master equation for the reduced density operator of the trapped atoms p is 
found by using standard techniques from the theory of open quantum systems [29] , and 
reads 

^ = ^ T i-i a \P a \ - r j _ 1 aja 1 p + H.c; (12) 

a similar equation is obtained for the hard-core bosons but replacing the boson operators 
by the corresponding spin operators. The effective interaction between different lattice 



sites is mediated by the Markov couplings (see |Appendix B[ ) 

r £ 



Ij_i = iexp (-ik L • rj_i) 



1 — erf I |j — 1| J — exp \—v-^ 



2X Q y V £ 



(13) 



where £ = l/d k with k = X _1 y 2|A|/w , and v = 1 (-i) for A < (A > 0). 
The error function is defined as erf (x) = (2/y/Tr) J Q X du exp (— u 2 ). This expression has 
been evaluated in the strong confinement regime coo — > oo and considering ki Xq 1 
(see Appendix B ). Note that the term '1 — erf (do |j — 1| /2X )' is basically zero for 



j 7^ 1, and therefore the £ parameter dictates the spatial range of the interactions as 



can be appreciated in figure |3j Finally, we remind that r = 4fi 2 y 2n\ A|/wg is the 
single-emitter decay rate. 



4-1. Extended Bose-Hubbard and spin models 

For k^ = and A < the Markov couplings are purely imaginary and have negative 
imaginary parts. Under this conditions, the master equation takes a Hamiltonian form 
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Figure 3. We plot the Markov couplings (real and imaginary parts at left and right, 
respectively) as a function of the distance between sites for A > 0. We have chosen 
do/Xo = 10 and = 0, and plotted 4 different values of £. It can be appreciated 
that this parameter controls the spatial range of the interactions. 



with effective Hamiltonians 

H e s = - ^2 h |r"a — il and H eH = - ^2 ft l r J~il °}°h ( 14 ) 

for non-interacting and hard-core bosons, respectively. The first Hamiltonian describes 
extended hopping in the lattice, a feature that has been recently shown to be helpful 
for achieving true incompressible Mott phases which otherwise can be blurred because 
of the additional slowly varying harmonic trap used to confine the atoms within the 
lattice [30] . The second Hamiltonian is equivalent to an extended ferromagnetic Ising- 
like Hamiltonian, and connects our system to extended spin models. Note that for j ^ 1 
the couplings have a Yukawa form 

l r j-i|^^77 ex pf-^V (15) 



|j - 1| /£ ' V £ 

and can even take a Coulomb form if £ is large enough. These kind of interactions 
(specially the Coulomb-like) are difficult to obtain with other techniques. 

These results show that in the trapping region the system under study could 
be useful for quantum simulation of condensed-matter phenomena requiring extended 
interactions. 

4-2. Dicke superradiance 

In this section we connect our system to superradiant Dicke-like phenomena both in 
atomic and harmonic oscillator samples (we will still take = 0). 

For A > 0, the Markov couplings are complex in general, and therefore the master 
equation of the system takes the form 
dp 1 



dt ih 



H d ,P 



V [p] , (16) 
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with a dissipation term given by 

V IP\ = ^2 7 J~' ( 2a iP a J ~~ ~ P°i<h) > ( 17 ) 
j,i 

having collective decay rates 

II 



7j _, = Re {T^} = r sinc j , (18) 

and a reversible term corresponding to inhomogeneous dephasing with Hamiltonian 

fl d = ^feA j _io}ai, (19) 



being 



erf ( — — li — II I — cos 



2X IJ 7 v £ 



IJ 



(20) 



The same holds for hard-core bosons but replacing the boson operators by the 

corresponding spin operators. 

For ( < 1 the Markov couplings do not connect different lattice sites, that is 
~ r 5j i, and the sites emit independently. On the other hand, when £ 3> M the 

collective decay rates become homogeneous, 7j_i ~ To, and we enter the Dicke regime. 

Hence, we expect to observe the superradiant phase-transition in our system by varying 

the parameter £. 

Note that in the Dicke regime the dephasing term cannot be neglected and connects 
the sites inhomogeneously with Aj^i ~ r £/ |j — 1|. This term appears in the optical case 
too, although it was inappropriately neglected in the original work by Dicke [10] when 
assuming the dipolar approximation in his initial Hamiltonian, and slightly changes his 
original predictions as pointed out in (3TJ [32] (see also [TT]). 

In the following we analyze the superradiant behavior of our system by studying 
the evolution of the total number of particles in the lattice = ^\ ^djaj^ and the 
rate of emitted atoms 

k 

in the last equality we have used that the total number operator byjbk + ^ dj dj is 
a constant of motion. 



4-2.1. Hard-core bosons: Atomic superradiance. Let us start by analyzing the case 
of a lattice in an initial Mott phase having one atom per site in the collisional blockade 
regime, which is the analog of an ensemble of excited atoms [9]. As explained in the 
Introduction, Dicke predicted that superradiance should appear in this system as an 
enhancement of the emission rate at early times [10], although this was later proved to 
happen only if the effective number of interacting emitters exceeds some threshold value 
[T2| [T3[ [H]: This is the superradiant phase transition. 
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Figure 4. Collective emission properties for an initial Mott state of hard-core bosons, 
(a) Time derivative of the rate at t = as a function of the range of the interactions 
£, see ([23]). The values M 3 = 8 (solid, blue), 27 (dashed, red), 64 (dashed-dotted, 
green), and 125 (dotted, yellow) are considered. It can be appreciated that there exists 
a critical value of £ above which the rate is enhanced at the initial times, (b) Rate as 
a function of time for M 3 = 27. The values £ = 0.01 (solid, blue), 0.5 (dashed, red), 1 
(dashed-dotted, green), and 10 (dotted, yellow) are considered. As expected from (a), 
the maximum of the rate is delayed above some critical £ value. Note that both the 
rate and its derivative have been normalized to the values expected for independent 
emitters, which are 4M 3 Fq and 2M 3 rg, respectively. 



In our system, the number of interacting spins is governed by the parameter £ (see 
figure IHJ), and the simplest way to show that the superradiant phase transition appears 
by varying it, is by evaluating the initial slope of the rate which can be written as[j| 



dt 



K 



-AM 3 T 2 



t=o 



1 



sine 



M 3 



(23) 



This expression has a very suggestive form: The term corresponding to the rate 
associated to independent emitters is balanced by a collective contribution arising from 
the interactions between them. In figure we show the dependence of this derivative 
with £ for various values of the number of sites M 3 . It can be appreciated that there 
exists a critical value of £ above which the sign of the derivative is reversed; hence, the 
rate increases at the initial time and we expect its maximum to be no longer at t = 0, 
which is a signature of superradiance. 

The time evolution of the rate for a cubic lattice with M 3 = 27 sites is shown 
in figure Hb for different values of £. We can appreciate how above some critical 
£ value the maximum rate of emission is delayed as expected. In order to find 
lZ(t) we have simulated the evolution equations for the coherences Cji = (cfjfri) and 

| Note that the evolution equation of the expectation value of any operator O can be written as 



m,l 



6, at. 



dt \ w / [dt 
and similarly the hard-core bosons in terms of the spin ladder operators. 



+ [&m,d])}, (22) 
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0.0 0.5 1.0 1.5 2.0 0.0 0.5 1.0 1.5 2.0 

T t T t 

Figure 5. Evolution of the total number of atoms in a lattice having M 3 = 27 sites 
for initial superfluid (a) and Mott (b) phases with N initial non-interacting atoms. 
The solid curves correspond to the limits £ -C 1 (grey) and £ 3> M (dark-blue) . Note 
how the 'evaporation time' is reduced for the initial superfluid state as £ increases (a). 
Equivalently, note how for an initial Mott state the atoms tend to stay in the lattice 
as £ increases (b). 



the populations Sj = (<t|), which we close by using the semiclassical approximation 
(4,ofoi) = (of) - 28 fl (olox); they read (E2} 

cj! = -4r Cji + ^2 r i-mCj m si + rj_ m c m iSj, (24a) 

m 

6i =-2j]r j _ 1 c jI + r?_ lClj . (246) 
1 

We have checked the validity of these semiclassical equations by comparing them with 
a direct simulation of the master equation for small number of sites in ID and 2D 
geometries; except for small quantitative deviations, they offer the same results. 

4-2.2. Non-interacting bosons: Harmonic oscillators superradiance. Let us analyze 
now the case of having non-interacting bosons in the lattice, which is equivalent to a 
collection of harmonic oscillators as discussed before. In previous works on superradiance 
this system was studied in parallel to its atomic counterpart [13], and here we show how 
our system offers a physical realization of it. We will show that superradiant effects can 
be observed in the evolution of the total number of atoms in the lattice, both for initial 
Mott and superfluid phases [§[ 

§ Let us note that a superfluid state with N excitations distributed over the entire lattice is more 
easily defined in the discrete Fourier-transform (DFT) basis 

j 

with q = (q x ,q y) qz) and q x ,y,z = 0, 1,2, ...M — 1, as the state having N excitations in the zero- 
momentum mode, that is, ISF)^ = (N\y 1/2 ft N |0). 



Quantum-optical phenomena with optical lattices 



14 



The evolution of the total number of atoms in the lattice is given by -see 

n T = -2^ 7j _iRe{(a}ai)}, (26) 
j,i 

and hence depends only on the real part of the Markov couplings. Therefore, we restrict 
our analysis to the dissipative term V [p] of the master equation (JTBjl . 

By diagonalizing the real, symmetric collective decay rates with an orthogonal 
matrix S such that S'pjTj-iS'qi = 7 P <5 pq , one can find a set of modes jc p = ^ S'pjOj j 
with definite decay properties. Then, it is completely straightforward to show that the 
total number of atoms can be written as a function of time as 

n T (t) = J] <c p c p (0)> exp (-2%t) • (27) 
p 

In general, 7j_i requires numerical diagonalization. However, in the limiting cases 
^ < 1 and £ ^> M, its spectrum becomes quite simple. Following the discussion after 
( T20l) . in the ^< 1 limit 7j_i = r <5ji is already diagonal and proportional to the identity. 
Hence, any orthogonal matrix S defines an equally suited set of modes all decaying with 
rate r . Therefore, if the initial number of atoms in the lattice is N, this will evolve as 

n T (t) = iVexp(-2r t), (28) 

irrespective of the particular initial state of the lattice (e.g., Mott or superfluid). The 
emission rate 1Z = 2T N exp (— 2r t), corresponds to the independent decay of the N 
atoms as expected in this regime having no interaction between the emitters. 

Let us consider now the opposite limit £ 3> M; in this case 7j_i = T V (j, 1), and 
the dissipative term can be written in terms of the symmetrical DFT mode only as 
V [p] = M 3 r (2/op/o - fofoP ~ io/op), see ((25). Hence, the DFT basis diagonalizes 
the problem, and shows that all the modes have zero decay rate except the symmetrical 
one, which has an enhanced rate proportional to the number of emitters. Therefore, 
starting with a superfluid state, the N initial atoms will decay exponentially with initial 
rate NM 3 T , that is 

n T (t) = iVexp (-2M 3 r t) . (29) 

On the other hand, if the initial state corresponds to a Mott phase, most of the atoms 
will remain in the lattice, as only the component which projects onto the symmetric 
mode will be emitted; concretely, from ( 1271) and ( 1251) the number of atoms in the lattice 
will evolve for this particular initial state as 

n T (t) =(n-^+^ exp (-2M 3 r t) . (30) 

Hence, according to this picture, superradiant collective effects can be observed in 
our system by two different means. Calling to the time needed to radiate the atoms in 
the absence of collective effects, one could start with a superfluid phase and measure 
this 'evaporation time' as a function of £; this should go from t for £ 1, to a much 
shorter time t /M 3 for £ ^> M (see figure Eh). Alternatively, one could start with a 



Quantum-optical phenomena with optical lattices 



15 



Mott phase, and measure the number of atoms left in the lattice in the steady state as 
a function of £; in this case, it should go from nT,steady = after a time t for £ <C 1, to 
n T ,steady = N- N/M 3 after a time t /M 3 for £ > M (see figure Eb). 

In order to find the evolution of nx we have simulated the equations satisfied by 
the coherences Cji = (ajai), which read -see (1221) - 

Cji = - ^ [lV m Cj m + rj_ m c ml ] . (31) 

m 

Note that in this case the equations are closed without the need of a semiclassical 
approximation, and hence they are exact. Note that they reproduce the analytic 
evolution of as given by (|28|) . ( 130|) . and ( |29l) in the corresponding limits (see figure 

Our results connect directly to those found by Agarwal some decades ago [T3] . 
Working with a Dicke-like model, he showed that if all the oscillators start in the same 
coherent state \a), the initial number of excitations, which in that case is given by 
iV = M 3 \a\ 2 , decays following ( J29l) . This is not a coincidence, but rather a consequence 
that, if N is large enough, a multi-coherent state of that kind is a good approximation 
of a superfluid state with that number of excitations. He also predicted that if the 
oscillators start in a number state, most of the excitations would remain in the steady 
state as follows from (|30|) . 

4-3. Directional superradiance 

The last phenomenon connected to quantum-optics that we want to discuss is how 
allowing 7^ the direction of emission of atoms can be controlled. 

This follows from the fact that by tuning the laser wave vector such that kt = ko 
while working in the A > regime, the master equation of our system (Tl2|) has exactly 
the same form as that of [18] (see also [H]), where the properties of light emitted by 
regular arrays of atoms were analyzed in terms of the rate and the direction of emission. 
Hence, we could expect the same phenomena to appear in our system, with the difference 
that in our case the emitted particles are atoms and not photons. Let us explain the 
results found in [18j and translate them to our system. 

We are particularly interested in the analysis performed in [18J of the emission 
properties of an initial symmetric Dicke state (a symmetrical spin wave in its notation). 
It was shown that the average number of photons emitted in the direction of the vector 
uq = (sin 9 cos 0, sin 9 sin 0, cos 9) , 9 and <fi being the polar and azimuthal spherical 
angles, respectively, can be written as 

1 (n) = E ex P P ^ " ■ r J-'] ( 32 ) 

where the total emission rate V is found from the normalization condition J dQI (Q) = 1 
(f2 is the solid angle). 

This expression was deduced in [18] in the M 3> 1 limit by diagonalizing the 
master equation with the use of periodic boundary conditions. The limit r c/L was 
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considered (c is the speed of light and L = Md is the length of the lattice in a given 
direction), that is, it was assumed that the emission time is larger than the time taken 
by a photon to leave the sample. Following the same steps, but taking into account the 
different dispersion relation of the emitted particles in our system, the same expression 
(1321) can be proved but now under the condition T <C v /L, where t>o = hk /m is the 
'resonant' speed of the emitted atoms. We will explain latter what this condition means 
in terms of the system parameters. 

In the optical case treated in [18], the angular distribution ( 1321) applies strictly only 
to the single excitation sector. In our case, this result is specially relevant for the case of 
non-interacting atoms trapped in the optical lattice in an initial superfluid state: since 
a superfluid state can be described as a state of N atoms in the completely symmetric 
state, see (l25|) . the results of [18] imply here that those N trapped atoms will emit iV 
free atoms with an angular distribution given by ( l32lh 

The sums in ( 132|) can be performed analytically, leading to the result 

sin 2 



m -*** n • V4>J (33) 



a=x,y,z 



where we have denoted by = k^/hi the unit vector in the 'laser' direction. This 
function has a diffraction maximum whenever the denominator vanishes, that is, for 
vectors m 6 1? such that u$] = k^ + 27r£m. In order to analyze directionality it is 
better to rewrite this condition as 7r£|m| = |k^ ■ m|, where rh = m/|m|. Then we see 
that for 7r£ > 1 only m = can satisfy this condition, and hence the atoms are mostly 
emitted in the laser direction. This is the directional regime. 

The emission rate T and the angular width of the atomic beam A8 must be 
evaluated numerically in general. However, approximating the diffraction peaks by 
Gaussian functions, a limit which is justified if the peaks are narrow enough (£/M <C 1), 
we can get the following expression for them [18J 

r = x r with X = n 3/2 M£ 2 , and A8 = £/M. (34) 

Hence we see that in addition to directionality, the rate is enhanced by a factor \- 

Let us finally explain what the initial approximation T ^ vq/L means in terms of 
the system parameters. Using (1341) . this can be written as 



^ 2 X /d 

which can be always satisfied by a proper tuning of the two-photon Rabi frequency. 

This analysis shows that our scheme can be used to observe directional, superradiant 
emission by starting with a superfluid phase and tuning the laser wave vector so that 
k L = k . 



Quantum-optical phenomena with optical lattices 



17 



5. Considerations about dimensionality 

So far we have allowed atoms in state |6) to move freely in 3D. However, the free motion 
of these atoms can be restricted to 2D or ID by coupling them to an additional ID or 
2D deep optical lattice, respectively. Under these conditions the momentum of the free 
atoms has components only along the non-trapped directions, say k 2 D = k x x + k y y and 
kiD — tx; from the point of view of the atoms in state \a), the dimensionality of the 
reservoir they are interacting through is therefore reduced to 2D or ID. Hence, only sites 
laying in the same z = const plane (for the 2D reservoir) or (y, z) = const tube (for 
the ID reservoir) would be connected through a two-point correlation function given by 
(ITT|) . but replacing G(t) by the correlation functions 

Q 2 Q 2 
C 2 d (t) = T j r exp (iAr) and G w (r) = —— exp (iAr) . (36) 

i 1 + 2 UJ o T ) (1 + ^qt) 7 

Gid coincides with the correlation function of atoms embedded in isotropic band-gap 
materials [2H [22] , while (?2D has no quantum-optical analog to our knowledge. 

Hence, with an additional trapping of the free atoms, our scheme could be useful 
to study low-dimensional reservoirs described by these correlation functions. 

6. Conclusions 

To conclude, let us briefly summarize what we have shown in this article. 

We have studied a model consisting of atoms trapped in an optical lattice (with 
trapping frequency u>o, lattice period do, and cubic geometry with M 3 sites) which are 
coherently driven to a non-trapped state through a Raman scheme having detuning A 
respect to the atomic transition, relative wave vector between the Raman lasers k^, 
and two-photon Rabi frequency Q. The limits of non-interacting bosons and hard- 
core bosons -where either one or zero atoms can be in a given lattice site- have been 
considered for the trapped atoms, while interactions involving the free atoms have been 
neglected. 

These are the main results that we have found (numbering follows the sections of 
the article): 

(121) We have first deduced the Hamiltonian of the system, showing that it is equivalent 
to that of a collection of harmonic oscillators or two-level systems for non- 
interacting atoms and hard-core bosons, respectively, interacting with a common 
radiation field consisting of massive particles. 

([3]) We have then studied the emission of free atoms when the lattice sites emit 
independently. We have found a phase transition at A = A — 4Q 2 /uq = from 
a regime (A < 0) where the trapped and free atoms form a bound state and 
there is a nonzero steady state population left in the lattice, to a radiative regime 
(A > 0) in which all the atoms leave the lattice eventually. We have identified the 
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Markovian and non-Markovian regimes of the emission, showing that the condition 
for Markovianity is |A| ^> 87tQ a /u)q. 

This behavior is analogous to the one predicted in [211 E2] for the spontaneous 
emission of atoms embedded in photonic band-gap materials. 

@ As for collective effects, we have studied them by using a reduced master equation 
for the lattice within the Born-Markov approximation, explaining under which 
conditions this may hold in our system. This master equation has shown that the 
number of interacting sites is controlled by a single parameter, namely £ = 1/g?o&o 



Now let us summarize the main collective effects that we have predicted to appear. 

(14.ip For A < and = 0, the evolution of the lattice follows an effective 
Hamiltonian: Extended hopping for non-interacting atoms and extended 
ferromagnetic interactions for hard-core bosons. The coupling between lattice 
sites has a Yukawa form, but approaches a Coulomb form as £ increases. 
This could be helpful for the simulation of condensed-matter phenomena 
requiring extended interactions. 

(14. 2p For A > and = 0, the master equation has a collective dissipation term 
in addition to the Hamiltonian part. In the ( < 1 limit this dissipation term 
corresponds to the lattice sites emitting independently with the same rate, 
while in the £ 3> M limit it takes a Dicke-like form, and therefore superradiant 
effects are expected to appear. Hence, we have argued that the superradiant 
phase-transition could be observed in our system by varying the £ parameter. 
Indeed, we have shown that both atomic and harmonic oscillator superradiance 
can be observed: 

( 14.2. ip When starting with a Mott phase having one atom per site in the collisional 
blockade regime (which is the equivalent of a collection of excited atoms), 
there exists a critical value of £ above which the emission rate is enhanced 
at the initial stage of the emission, and its maximum occurs for t ^ 0. 
This delay of the maximum rate was suggested by Dicke in his seminal work 
[TO] , and was shown to be a signature of the superradiant phase-transition 
for atomic samples [121 H31 EH] • 

(I4.2.2p If the lattice atoms are non-interacting, there are two equivalent ways of 
observing superradiant behavior. One can start with a superffuid phase, 
and measure the time needed to radiate the lattice atoms, which should 
be reduced by a factor M 3 when going from the ^ C 1 to the £ ^> M 
limit. Alternatively, one can start with a Mott phase having N atoms in 
the lattice, and measure the number of trapped atoms left in the steady 
state; this should go from zero in the £ < 1 limit, to N - N/M 3 in the 
£ > M limit. 

This behavior is consistent with the results found in [13] for the 
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spontaneous emission of a collection of harmonic oscillators, and hence 
our scheme provides a physical realization of this system. 

( 14. 3 p The last phenomenon linked to quantum-optics that we have shown is that 
superradiant, directional emission can be obtained in our system. In particular, 
we have proved that by starting from a superfluid phase distributed over enough 
lattice sites (M 3> 1), matching kt to ko, and working in the 7r£ > 1 regime, the 
atoms are emitted mostly within the direction. If the condition £/M <C 1 
is also satisfied, the beam has been shown to have angular width AO — £/M 
and a rate enhanced by a factor \ = 7r 3 / 2 M£ 2 . 

These results are in correspondence with those of [T7J [18] , where the emission 
of photons by regular arrays of atoms was studied. 

05]) We have finally discussed what would change if the motion of the free atoms is 
restricted to 2D or ID by using an additional optical lattice. We have shown 
that in this conditions the correlation function of the reservoir through which the 
lattice sites are interacting takes the form G(r) oc (1 + iojqt /2) d l 2 , where d is the 
dimensionality {d — 1 and 2 for ID and 2D, respectively). For ID this corresponds 
to the correlation function of the radiation field in an isotropic photonic band-gap 
material [21, 22J, while for 2D it has no quantum-optical counterpart. 

All the proposed scenarios can be reached with currently available technology, and 
hence we believe that this work provides a solid base for considering optical lattices as 
quantum simulators of quantum-optical systems. 
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Appendix A. Extracting information from the non-Markovian equation (J7j) 

In this Appendix we explain how to manipulate ([Z]) by using Laplace transform 
techniques in order to understand the emission properties of the system in the limit 
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of independent lattice sites. 

Making the Laplace transform of (JTJ) the solution for the amplitude A in Laplace 
space is found as [jj 

1 (A.3) 



A(s) 
V ' s + G{s) 

On the other hand, the Laplace transform of the correlation function (jSJ) is 

An 2 



G(s) 



-i + e 



_2(A+i£) /27r(A + is) 



CJ 



1 + erf U 



'2(A + is 




or to first order in (A + is) /ojq 



G r 



s = -- 



UJ 



+ a(l + i) V27r(s-iA). 



(A.4) 



(A.5) 



Note that as we assume that s/uq <C 1, this transform cannot describe properly times 
below oJq , and this is why sometimes this limit receives the name 'strong confinement', 
as it needs Uq — > oo or otherwise the initial steps of the evolution are lost in the 
description. Unlike in the general case where the transform of the correlation function 
has a really complicated form (IA.40 . in this limit it is possible to make the inverse Laplace 
transform of A (s). Hence, in the following we consider the two cases separately, as the 
first case requires special attention in order to extract results. 



Appendix A.l. Evolution of the population in the strong confinement limit 
The Laplace transform of A can be inverted as [33] 

i re+ioo -r re+ioo st 

e - iAt A (t) = — / ds e st A (s + iA) = — / ds . , (A.6) 



^" 1 J t— ioo 

z,nLj e _ ioo s + iA + a (1 + i) V2ns 
where we have used the 'strong confinement' form of the correlation transform (|A.5[) . 
and e > must be chosen so that all the poles of A (s + iA) stay to the left of the 
Re {s} = e line. The multivaluated character of y/s forces us to define a branch cut in 
the complex-s space; we chose to remove the branch arg{s} = it, hence defining the 
phase of s in the domain ]— it, 7r[. Then, choosing the integration contour as that in 
figure I AH we can write the needed integral as 



I dsK(s) = 2iriJ2 R i- lim / dsK (s) + 



dsK(s) 



(A.7) 



R-too 



We define the Laplace transform of a function / (t) as [33] 



f(s) = £ s [/(*)] = 
which satisfies the properties [33] 



dte~ st 



f{t) =sf(s)-f(0) and C s 



/(*), 



dt'g(t-t')f(t') 



(A.l) 



(A.2) 




Figure Al. Integration circuit in the complex-s plane for the evaluation of the inverse 
Laplace transform. 

being Rj the residues of the kernel K (s) = A (s + iA) exp (st) at its poles [33]. The 
integrals over the paths C u , Cd, and C c are zero when r — > and R — > oo. 

Making the variable change s = e m x and s = e~ m x in the integrals along 1Z U and 
IZd, respectively, we get 

r+oo /— -xt 

A^ae~ i7T/4 I dx ^ , (A.8) 

{^—x + iA j + i4na 2 x 



lim 

R— too 



dsK(s) + / dsK(s) 
Un u Jn d 



o 



which is not analytical, but can be performed numerically very efficiently 

On the other hand, the poles of K (s) are easily evaluated by making the change 
s = ix 2 , hence turning its denominator into x 2 + l^pnax + A whose roots are 



x± = - v^Fa ± V na 2 - A = \x± | e iv± ; (A.9) 

therefore the poles are simple and read s± = ix 2 ± . In order to evaluate the residues, it 
is convenient to write the Kernel of (1A.6jl as 

g w=n f_ + e,C> ^ (Aio) 

j=± j 
from which the residues are easily found as 

R + = \± -± e s+ * an d R_ = ^— — — e s "*. A.ll 

v /s^ - e 17r /4 x _ ^sZ - e m /*x + v ' 

The problem now is that the phase of s± must be defined in the interval ]— 7r,7r[ 
(remember the branch cut in figure IA10 . and hence we have to analyze their phases 
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as a function of the phases of x±, which we take in the interval [0, 2ir] as a convention. 
In terms of the parameters A and a, we can define 3 regions: 



In the region A < we have <p + = and ip_ = ir. Hence, we have to write 

s+ = e [7r/2 \x+\ 2 and s_ = e i7r/2 |x_| 2 e 2i " x e~ 2i ", (A.12) 
that is, we have to decrease the phase of s_ by 2tt to keep it within the interval 



-7r, 7r[; with s± written in this way, the residues read 



2 


\x+ 




\ x +\ 


+ | 





Axi_t 



and R_ = 0. 



(A.13) 



In the region < A < na 2 , we have ip± = tt, and hence, we have to decrease the 
phase of both s± by 2n to keep them within the interval ]— 7T, 7r[; therefore both 
residues are zero R± = 0. 

In the region A > na 2 , we have <p + G [tt/2, ty], and = [7r, 37r/2]. Hence, we need 
to write 



e in/2 \x±\ 2 e 2iv± xe- 2in and 



xe" 



-4i7T 



(A.14) 



that is, we decrease the phase of s + by 2n and the phase of s_ by Att, so that they 
both stay in the interval ]—tt, 7r[, and obtain 

2x_ 



and i?_ 



(A.15) 



Introducing these results for the residues and the integral (1A.8I) in (1 A. 71) . we arrive 
to the expression that we wanted to prove for the population A (t), see 0. 



Appendix A. 2. Steady state population in the general case 

From the last discussion it should be clear that the long term population is dictated 
solely by the residues of K (s) = A (s + iA) exp (st). Moreover, only the residues 
corresponding to pure imaginary poles will remain in the t — > oo limit, and hence 
the steady state population will be given by 

2 



\A(t -> oo) 



E 



R i 



(A.16) 



Sj/Re{sj}=0 

Then, in order to know the steady state population in the finite trap case, we just need 
to find the pure imaginary poles of K (s) with the transform of the correlation function 
given by ( ]A.4j) . By writing s = ix with x6l, this amounts for solving the equation 



x + A + 4fi 2 exp ( — ) 

w y 



2nx 



orf ! \ I — 



0. 



(A.17) 
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This equation cannot be solved analytically, but it can be efficiently done 
numerically. On the other hand, once we know the solutions Xj of this equation, we 
can use them to evaluate their associated residues as [33] 

-l 



R 3 



'ix ,• t 



ds 



s + iA + G (s + iA) 



(A.18) 



This is the method that we have used to compute the results shown in figure EH. 



Appendix B. Evaluation of the Markov couplings ( 113ft 

In this appendix we show that in the strong confinement regime the Markov couplings 
have the form ({TBI . After the Born-Markov approximation is carried, the effective 
coupling between different sites has the form 

9k ex P [i ( k - k • rj-i] 

£->0 + 



rj 1= / drG j _ 1 (r)= hm V 



h 2 



e + iA h 



(B.l) 



In the limit Ul <C Xq 1 (in which we can remove k^ from -see (jlj)-), taking the 
continuous limit for the momentum, and performing the angular part of the integral 
over momenta by assuming that rj_i is within the z-axis, this expression can be reduced 
to 



r.i-i 



8x n 2 



,- ik i r j-i 



lim ( A + ie 

e->0+ 



d/csinc (krj-i) exp (— X$k 2 ) 

sine (A;rj_i) exp (— Xq k 2 ) 



(B.2) 



The first integral is easily evaluated as 

dk sine (kr- } -\) exp (— Xq k 2 ) 



71 



2r 



-erf 



j i 



A + ie - 

2m 



r i~l 



2X r 



(B.3) 



The second integral is not analytical in general. Nevertheless, in the strong confinement 
regime both |A|/o;o and Xq/oIq are small, and hence the Gaussian function in the Kernel 
can be approximated by one, so that the remaining integral reads 



dA;^ 



sine (kr 



j-U 



7T 



A + i£-§± 2r j _ 1 (A + ie) 



1 — exp i 



r Hi /2(A + ie) 
X 



(B.4) 



Introducing these integrals in (1B.2[) and talking the e — > limit we arrive to expression 
( TT3|) for the Markov couplings. 

Note that e is more than just a mathematical artefact for regularizing the integral. It accounts for 
the decoherence channels not taken into account in the model, which make the correlation between 
lattice sites decay as Gj_\(t) exp(— er). Even if these channels appear at a time scale larger than the 
relevant times (e <C 1), they remove the zeros of the kernel's denominator, showing that its singularities 
are not physical. 
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